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Summary 

Chaos  in  Classical  Nonlinear  Fields 


We  studied  chaotic  behavior  in  a  high-dimensional  periodic  lattice  version  of  the 
classical  Hamiltonian  phi-four  field  theory.  Both  single  and  double  well  potentials  were 
considered.  We  examined  the  existence  of  a  global  stochasticity  threshold  that  varied 
with  energy  and  initial  conditions.  The  model  was  discretized  using  an  algorithm  due 
to  Hirota  that  guarantees  stability  and  energy  conservation.  Chaotic  behavior  was 
diagnosed  using  the  Lyapunov  exponent,  with  additional  information  from  space-time 
profiles,  Fourier  pow'er  spectra,  and  phase  space  plots.  We  found  long  time  scales  which 
made  it  difficult  to  distinguish  between  asymptotically  chaotic  and  integrable  behavior. 


1.  Background  and  Motivation 


The  importance  of  understanding  infinite  (or  high)  dimensional  nonlinear  systems 
is  clear:  we  are  surrounded  by  fluids,  plasmas,  solids  and  other  continuous  (or  very  high 
dimensional)  systems  whose  behavior  is  critical  to  us.  In  some  cases,  when  dissipation 
is  present,  the  long-time  behavior  of  such  systems  may  be  adequately  described  by  low 
dimensional  models.  But  for  Hamiltonian  systems  with  no  asymptotic  description  on 
attractors  possible,  and  where  regular  and  chaotic  regions  can  be  strongly  intermixed, 
we  usually  must  seek  understanding  and  intuition  in  a  different  way. 

Another  reason  for  the  importance  of  high  dimensional  Hamiltonian  systems  is 
their  relevance  to  basic  questions  in  statistical  mechanics,  including  the  nature  of  the 
approach  to  equilibrium  and  the  equipartition  of  energy.  Recent  equipartition  studies1 
in  large  Hamiltonian  systems  containing  no  a  priori  natural  time  scale,  suggest  that 
very  long  times  are  necessary  for  some  modes  to  reach  equilibrium. 

In  the  present  work  we  studied  a  non-integrable  classical  Hamiltonian  field  theory, 
the  <f> 4  model  (which  will  be  described  below).  This  arises  physically  as  a  useful  simple 
model  in  condensed  matter  systems  (e.g.,  polyacetylene,  ferroelectrics)  and  elementary 
particle  physics.  But  its  basic  significance,  for  our  purposes,  is  as  an  example  of  a 
nearly-integrable,  high-dimensional,  nonlinear  Hamiltonian  system. 

Our  working  hypothesis  was  that  chaotic  behavior  is  to  be  expected  in  a  non- 
integrable  system  for  some  values  of  the  parameters  and  types  of  initial  conditions. 
We  expect  that  at  sufficiently  high  energy  (or  energy  density)  global  stochasticity  will 
occur  in  both  the  single-  and  double-well  potentials,  with  the  two  cases  being  rather 
similar  in  the  high  energy  regime.  Since  only  the  double  well  allows  solitary  wave  (and 
near-breather)  solutions,  we  expected  the  low  energy  behavior  of  the  two  models  to  be 
rather  different. 

Previous  numerical  studies'  have  found  conflicting  results,  some  finding  no  chaos 
or  chaotic  behavior  which  appears  as  the  energy  (or  energy  density)  increases  and  then 
disappears  at  higher  energies  How  the  partial  differential  equation  is  discretized,  what 
initial  conditions  are  considered,  choice  of  lattice  parameters,  and  especially  the  length 
of  integration  times  are  expected  to  be  crucial  to  the  results. 


I 


2.  The  (j) 4  Model 

The  Lagrangian  and  Hamiltonian  densities  of  the  (one  space  dimensional)  field 
theory  are: 


where 


The  equation  of  motion  for  this  potential  is 


d7(f>  d7<t> 

dt2  dx 2 


+  “f-  \<j>  —  o 


This  model  can  be  regarded  as  a  truncation  of  the  exactly-integrable  sine-Gordon  (/x  < 
0)  and  sinh-Gordon  (/x  >  0)  equations.  With  /x  <  0,  V(<f> )  is  a  double- well  potential 
(leading  to  spontaneous  symmetry  breaking  in  the  corresponding  quantum  field  theory); 
the  classical  theory  then  allows  solitary  wave  solutions. 

It  is  possible  to  remove  the  (f>2  and  <j>A  coefficients  in  F(<£)  by  scaling  x,t,  and  4>  : 


t  - 1 h_1/2< 


(^y 


Thus  the  equation  of  motion  can  always  be  written  in  the  form 

4>tt  -  ±  <t>  +  <f>3  =  0  (<f>t  =  ,  etc.)  [1] 

In  this  form  we  see  that  the  only  parameters  which  can  be  varied  are  the  total  energy 
and  sign  of  the  quadratic  term  in  the  potential  (corresponding  to  the  single  and  double 
well  cases).  For  the  lattice  version  of  the  model  there  are  additional  parameters,  the 
number  of  lattice  points  and  length  of  the  lattice. 
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The  <t>4  model  has  two  known  conserved  quantities: 


E  = 


P  =  - 


J  4>x<t>tdx 


The  absence  of  an  infinite  number  of  conserved  quantities3  is  also  supported  by  studies 

of  (f> 4  solitary  wave  scattering  which  show  that  the  initial  kink  and  antikink  shapes  need 
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not  survive  the  collision. 


3.  Progress 

This  work  began  as  a  Master’s  thesis  of  D.  Bradley-Hutchison,5  supervised  by  one 
of  us  (HKS).  He  used  a  300  site  spatial  lattice  (with  periodic  boundary  conditions), 
a  lowest  order  finite  difference  approximation  to  the  equation  of  motion,  and  initial 
conditions  consisting  of  a  sharply  peaked  static  gaussian  at  the  center  of  the  lattice 
A  simple  diagnostic  was  used  to  measure  the  separation  of  neighboring  trajectories  in 
phase  space: 


m)}7  = 


i 

L 


/  dx  {  4>{t)  -  <£(t)]  +  [<£t(0  - 


where  (f>  and  (f)  are  initially  infinitesimally  close  trajectories  .  For  an  integrable 
system  we  expect  D(i)  ~  t  for  large  t,  while  for  a  chaotic  system  we  expect  In  D(t )  ~~  t 
asymptotically. 

Bradley-Hutchison  studied  the  separation  d’stance,  D(t),  as  the  total  energy  was 
varied  over  the  range  43  <  E  <  1865.  He  found  that  D{i )  grows  linearly  for  “low”  ener 
gies  and  “short”  times,  corresponding  to  typical  integrable  behavior,  and  that  lnD(/) 
grows  linearly  for  “high”  energies,  corresponding  to  chaotic  behavior.  For  “moder¬ 
ate”  energies,  there  can  be  a  sudden  onset  of  chaotic  behavior  after  “long”  times.  In 

the  strongly  chaotic  cases,  “saturation”  is  generally  observed;  i.e.,  In  D(i)  reaches  a 
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maximum  value  at  some  t  value  and  remains  constant  thereafter  (see  Fig.  1). 
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These  preliminary  results  were  interesting  enough  to  warrant  a  more  detailed  and 
sophisticated  study.  There  are  two  primary  components  in  the  calculation:  the  dis¬ 
cretization  algorithms  for  the  partial  differential  equation,  eq.[l],  and  the  selection  of 
diagnostics  for  the  spatial-temporal  behavior.  Ideally  we  would  like  to  discretize  in  such 
a  way  that  any  known  conserved  quantities  for  the  continuous  system  are  preserved, 
and  also  in  such  a  way  that  stability  is  guaranteed  as  the  spatial  and  time  step  sizes 
are  varied  in  the  integration  algorithms. 
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We  chose  a  discretization  scheme  due  to  Hirota  that  guarantees  stability  and 
conservation  of  energy.  The  discretized  Hamiltonian  density  is 

n  =  \e2  +  j(A2,*)2  +  V(4). 

The  corresponding  first  order  equations  of  motion  are 

a  t4>  =  n  tp 

Atp  = 

where  the  central  difference  operator,  A,  and  averaging  operator,  II,  are  defined  as 

a  ,/(*)  =  ![/(*  + 5) -/(*_!)] 
n./(*)si[/(*  +  ^)  +  /(*-|)] 

with 

S2zf{x)  =  AtUxf{x) ,  etc. 

Combining  the  first  order  equations,  we  obtain 

a  ict>  =  Alxn2t<j>  -  nt& ,v(*)| 

The  implicit  nature  of  the  Hirota  integration  algorithm  means  that  in  genera]  we 
must  solve  a  set  of  coupled  nonlinear  equations.  When  the  equations  are  linear,  the 
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resulting  equations  can  be  solved  using  techniques  for  inverting  a  triuiagonal  matrix 
(modified  for  the  periodic  boundary  conditions  that  we  are  imposing).  With  the  non¬ 
linear  terms  included  in  the  evolution  equations,  we  use  an  explicit  scheme  as  an  ap¬ 
proximation  to  the  nonlinear  term  and  then  iterate  the  inversion  until  the  field  becomes 
self-consistent. 

In  our  work  we  used  a  simple  measure  of  chaotic  behavior,  the  maximum  Lyapunov 
exponent,  to  search  for  the  onset  of  global  stochasticity.  This  exponent  is  a  function  of 
the  linearized  equations  of  motion  which  are  in  turn  dependent  on  the  full  equations 
of  motion.  On  a  discrete  lattice,  the  equations  for  the  field  become  a  set  of  N  coupled 
oscillators  governed  by  the  equations  of  motion  of  the  form 

d(f 
It 

where  <f>  is  a  2 N  dimensional  vector,  N  components  are  the  field  and  the  other  N  are 
its  time  derivative.  We  linearize  the  equations  of  motion  about  our  initial  trajectory 
by  introducing  the  variable  0. 

r  a  _  ,• _  0(00  +  *0,0  ~  0(00,0 

0(00,  i)  =  »m  - TTr, - 

\6<j>\ 


W) 


where  4>o  indicates  the  initial  values  for  the  variables.  These  new  variables  evolve  in 
time  according  to  the  the  linear  equations 


c/0, 

dt 


;  94> j 


0, 


0(00. <) 


These  time  evolution  equa'-  :;s  f  t  0,  are  a  set  of  2N  coupled  ordinary  differential 
equations  and  were  integral  i  ■•i-ing  a  fourth  order  Runge-Kutta  algorithm.  Finally, 
the  maximum  Lyapunov  ex;  -  nt  is  given  by 


lim  -  In 

1— >oo  t 


110(^0,011 

l|0(0O,O)|| 


If  there  is  sensitivity  to  initial  •  millions,  nearby  trajectories  diverge  exponentially  and 
the  Lyapunov  exponent  is  positive.  For  integrable  systems  the  divergence  of  nearby 
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trajectories  is  linear  and  the  finite  time  Lyapunov  exponent  decreases  toward  zero  as 
t-1.  D(t)  and  the  Lyapunov  exponent  measure  essentially  the  same  property  of  the 
system;  however  the  Lyapunov  exponent  procedure  is  preferable  in  that  it  takes  the 
limit  of  the  initial  separation  going  to  zero. 

Figures  2-4  show  some  of  our  results.  In  these  calculations  (performed  on  the  UNH 
VAX  8820  and  the  San  Diego  SDSC  CRAY  X-MP)  we  imposed  periodic  boundary 
conditions  on  discrete  lattices  with  128  or  256  points.  Energy  was  conserved  to  6 
or  7  digits  at  each  time  step.  The  initial  conditions  examined  to  date  include  static 
gaussians  of  varying  amplitude  centered  on  the  lattice  and  pure  sinusoidal  modes. 

As  a  general  summary  of  these  runs,  we  note  (for  initial  static  gaussians):  (a)  there 
is  little  dependence  of  the  Lyapunov  exponent  on  the  number  of  lattice  points  (for  N 
above  100);  (b)  for  energies  up  to  80  in  both  the  single  and  double  well,  we  see  the 
Lyapunov  exponent  decreasing  as  t_1,  indicating  non-chaotic  behavior;  (c)  at  E  ~  300 
the  Lyapunov  exponent  approaches  a  positive  constant  value,  indicating  chaos;  here 
the  double  well  has  a  larger  exponent  than  the  single  well;  (d)  for  all  energies  between 
E  =  300  and  E  =  2600  the  asymptotic  t  behavior  is  not  yet  clear,  since  the  data  shows 
regions  of  flatness  as  well  as  sharp  decreases,  (c)  for  high  energies  (2600  and  17356)  the 
system  again  has  a  decreasing  Lyapunov  exponent.  Thus,  from  our  search  in  energy 
space  with  static  gaussian  initial  conditions,  it  appears  that  there  might  be  an  onset 
and  later  disappearance  of  global  stochasticity.  The  study  of  this  question  became  the 
focus  of  our  work  at  the  end  of  the  grant  period. 

Code  development  for  all  the  additional  useful  diagnostics:  space-time  profiles, 
power  spectra,  and  phase  space  plots  has  also  been  completed. 

In  examining  the  Lyapunov  exponents  over  very  long  times,  we  became  aware  of 
the  sensitivity  of  our  results  to  the  details  of  the  discretization  procedure.  For  this 
reason  it  became  difficult  to  be  certain  if  the  Lyapunov  exponent  differed  significantly 
from  zero  and  thus  whether  the  system  was  chaotic  or  not.  We  are  currently  devising 
more  sophisticated  statistical  tests  to  be  applied  to  our  Lyapunov  exponent  data  in 
order  to  distinguish  between  zero  and  non-zero  values. 
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4.  Communications  and  Publications 

This  grant  has  resulted  in  the  publication  A  Numerical  Study  of  Chaos  vn  the  One- 
Dimensional  d>*  Model,  Douglas  A.  Bradlev-Hutchison  and  Harvey  K.  Shepard,  Physica 
Scnpta  40,  731  (1989),  and  the  presentation  of  a  poster  at  Dynamics  Days  Conference, 
Houston,  Texas  in  January  1989. 
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lyopunov  exponent 


Varying  Lottice  Size 
nx=64,  128,  256 
initiol  gaussion:  A==9,  W=50 
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Single  Well;  W=50;  nx=128 


